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Abstract 

A simple lattice model is used to study compaction in granular media. As 
in real experiments, we consider a series of taps separated by large enough 
waiting times. The relaxation of the density exhibits the characteristic inverse 
logarithmic law. Moreover, we have been able to identify analytically the rel- 
evant time scale, leading to a relaxation law independent of the specific values 
of the parameters. Also, an expression for the asymptotic density reached in 
the compaction process has been derived. The theoretical predictions agree 
fairly well with the results from the Monte Carlo simulation. 

PACS numbers: 81.05.Rm,05.50.+q,81.20.Ev 
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I. INTRODUCTION 



One of the characteristic complex behaviors exhibited by granular materials is com- 
paction [l|-|4| . It can be roughly defined as the density relaxation of a loosely packed system 
of many grains under mechanical tapping or vibration. Granular compaction is important 
to many industrial applications related with the production and manipulation of a wide 
variety of systems composed by many macroscopic particles or grains [Q. In the last few 
years, a series of experiments have been carried out trying to identify the physical principles 
underlying granular compaction fHHH . Starting from a loosely packed initial configuration, 
systems of monodisperse glass beads were tapped vertically. The waiting time between suc- 
cessive taps was large enough to allow the system to relax, so that the beads were at rest 
before the next tap started. The time evolution of the density towards a steady state has 
been analyzed, and it has been shown that it can be accurately described by an inverse 
logarithmic law with four adjustable parameters, whose values depend only on the tapping 
strength measured by the peak acceleration of a tap. The logarithmic relaxation has been 
found in many different models ||-[|] suggesting that such a behavior is quite general 0. 
Although several mechanisms have been proposed to explain the behavior observed in the 
experiments, a fully satisfactory theory is still lacking. 

Here we consider a one-dimensional model simple enough as to allow some detailed 
calculations. One of our aims was to try to identify the relevant time scale over which 
the relaxation (compaction) of the system takes place. This is the first step in the search of 
general laws governing the physics of densification. A main difficulty in studying compaction 
is that there are two different series of elementary processes involved in the experiment. The 
system is submitted to taps or pulses separated by time intervals for which the system is 
allowed to relax freely. The initial state for each tap is the final state from the previous 
relaxation. Both processes, tapping and free evolution, must be considered in detail, and 
they are rather different from a physical point of view. For instance, while the duration of 
the pulse is clearly a relevant parameter of the problem, the free relaxation is assumed to 
last by definition until the system gets trapped and it is at rest. Quite interestingly, the 
experiments have shown that it is useful to measure time by the number of pulses applied to 
the system. Another central question is how much settling will occur for a given vibration 
intensity, and also if the stationary value of the density depends on the initial configuration. 
How these facts appear in our model and which is the role played by the duration and 
amplitude of the pulses are points we will address here. 

The plan of the paper is as follows. In the next Section the model will be presented. It 
consists of a lattice whose sites can be occupied by particles. The dynamics is formulated by 
means of a master equation and it is facilitated, in the sense that the rates of adsorption and 
desorption of a particle are proportional to the number of particles in the nearest neighbor 
sites. The model can be exactly solved in the no desorption limit, which corresponds to the 
very low temperature limit. The solution is obtained in Sec. |T|and describes the evolution of 
the system without external perturbation. Therefore, it will be used to study the relaxation 
of the system towards a metastable state between pulses. In spite of the simplicity of the 
system, the general solution for arbitrary strength of the external energy source is rather 
complicated. We have considered the limit of short duration of the taps, not only because 
of mathematical convenience, but also because it seems to be the limit in which the time 
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scales involved in the problem become well separated. 

The sequence of taps and free relaxation processes, i. e. compaction is the subject of 
Sec. [TV]. An expression for the density after the n + 1 tap in terms of the density and the 
probability distribution of two holes separated by a site after the previous tap is derived. 
Although this relation does not provide an explicit expression for the evolution of the density, 
it allows to identify the relevant time scale, which turns out to be proportional to the duration 
of a tap times a parameter measuring their strength. Curves describing the density evolution 
of systems starting from the same initial state but corresponding to different values of the 
parameters are shown to be the same when plotted as functions of the scaled time. Besides, 
the single scaled curve is very well fitted by the inverse logarithmic law known from real 
experiments. 

In the limit of many taps, the density reaches a steady value that is discussed in Sec. 
[V[ By using a pair approximation it is found that the steady density is proportional to the 
time relevant parameter mentioned above. This prediction agrees well with the numerical 
results from the simulations. To the best of our knowledge this is the first time an explicit 
expression for the asymptotic density limit in a discrete tapping process has been analytically 
derived from a model for compaction. Finally, the last section contains some final remarks 
and comments. 



II. DESCRIPTION OF THE MODEL 

We consider a one-dimensional lattice with N sites. Each site can be either occupied by 
a particle or empty. A configuration of the system is specified, for instance, by giving an 
ordered sequence of N particles and holes. Let us introduce a set of variables m = {m^; i = 
1, 2, • • • , N}, such that m>i vanishes if there is a particle at site i, while it takes the value 1 
if there is a hole, i. e. the site i is empty. 

The dynamics of the system is defined as a Markov process and formulated by means of 
the master equation for the conditional probability pi/i(m, t\m', t') of finding the system in 



the configuration m at time t, given it was in the configuration m' at time t' <t fllP 




Of 



Pi /1 (m,t|m',t / ) = [w i (R i m)p 1/1 (R i m,t\m',t') - W r i (m)pi /1 (m,t|m',t / ) 



where = {mi, • ■ ■ , Rirrii, • ■ ■ , m^} with Rivrti = 1 — rrii, i. e. R^m is the configuration 
obtained from m by changing the state of hole or particle of site i. The above equation is 
to be solved with the initial condition 

JV 

pi/i(m,t'|m / ,t / ) = 6 m>m ' = n<Wm<- ( 2 ) 

8=1 

The one-time distribution 

p(m,t) = ^pi/i(m,i|m',0)p(m',0) (3) 

m' 

also obeys Eq. (ffl), although now the initial condition must be given in each specific situation. 
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The possible elementary processes occurring in the system are the adsorption of a particle 
on an empty site from a surrounding bulk and the desorption of a particle from the lattice to 
the bulk. The probability that an adsorption attempt be made on site i in the infinitesimal 
time interval between t and t + dt is k + dt. Of course, a particle can be adsorbed only if the 
site is empty In the same way, the probability per unit of time that a given particle try to 
leave the lattice is Both processes are restricted in the following way. A particle can 
be adsorbed on or desorbed from a site only if at least one of its nearest neighbor sites is 
empty. More precisely, the probability rate for the events is proportional to the number of 
nearest neighbor holes. This condition tries to model naively the short-ranged geometrical 
constraints that make structural rearrangements difficult in a granular material. Therefore, 
we assume that the transition rates are given by 

Wi{mS> ) m i = l) = k + , Wi(m w , = 0) = , (4) 

with = {mi, • ■ • , mj_i, mj+i, • • • , tun}. A similar kind of facilitated dynamics has been 



used previously in the context of Ising models ||11|| . Although we restrict ourselves here 
to the one-dimensional case, the model can be formulated for arbitrary dimension. Let us 
introduce new constant parameters v and e by 

" = *- + *+' e = FTV (5) 

The constant v has the dimensions of a frequency and e is a dimensionless parameter defined 
in the interval < e < 1. For e = 1 no particle is adsorbed by the system, while for e = 
desorption processes do not occur. In terms of these parameters, Eqs. (H) can be written 
together as 

Wi(m) = |(mi_i + m i+1 )[e + m^l - 2e)]. (6) 
The ratio of the desorption transition rates to the adsorption ones is 

Wi(m<- i \m i = 0) k_ e 



Wiiva^.rrii = 1) k+ 1 - e 



(7) 



where the last equality defines the parameter x. 

The stochastic process we have formulated has a steady one-time distribution of the form 



i 



and the density of holes (average number of holes divided by the total number of sites N) 
in the steady state is 



(mi) st = J2 m iP8t(m) = — — = e, (9) 

1 + x 



and, consequently, the steady density of particles is p st = 1 — e. As x increases, the equilib- 
rium density of particles decreases. 
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The quantity x(e) can be related with a temperature parameter T by defining an energy 
E(m) for the system. A possible choice is 



N 



E{m) =e J2m, (10) 



where eo is a constant fixing the energy scale. If now the distribution given by Eq. (]§]) is 
identified with the equilibrium canonical distribution, it is easily obtained that 

x = e"^ , (11) 

with (3 = (ksT)^ 1 , ks being the Boltzmann constant. Thus the limit x — > oo (e — > 1) is 
equivalent to T — > 0~ and the limit x — ■> (e — ► 0) to T — > + . A purely random distribution 
of particles and holes corresponds formally to the equilibrium distribution for e = 1/2 or 
T — > oo. 

The transition rates given in Eq. (|7|) define an irreducible Markov process for e > 0, 
except for the state with all the sites occupied by particles that can not evolve. In the 
limit N — > oo, the probability of this state is negligible, and all the solutions of the master 



equation relax to the steady distribution given by Eq. ||10|| . The situation is different in 
the no desorption limit e = 0. The density of particles cannot decrease, and all the states 
of the system having every hole surrounded by two particles are absorbent; no evolution is 
possible from them. 

Our one- dimensional lattice model can be regarded as a very simple picture of an hori- 
zontal section of a real granular system, near the bottom of the container. Consider first the 
freely evolving case. In a real granular medium, particles cannot go up due to gravity. They 
can only go down, as long as there is enough empty space in their surroundings. Therefore, 
the packing fraction grows until the hard-core interaction prevents more movements of par- 
ticles, and a mechanically stable configuration is found. This situation is naively resembled 
by the evolution of our model in the no desorption limit, e = 0. Starting from a given 
configuration of particles and holes, the system evolves by means of adsorption processes, 
occurring on those sites having at least one nearest neighbor hole. This leads to a monotonic 
increase of the density until all the holes become isolated, i. e. surrounded by two particles. 

Next, suppose a granular system submitted to vertical vibration. During the vibration, 
particles belonging to a low horizontal section can go up, making the local packing fraction 
decrease. The hard-core repulsion is also fundamental in the vibrated case, since particles 
always need enough free volume close to them in order to move. In our model, these pulses 
are introduced by allowing particles to be desorbed, but the dynamics is "facilitated". A 
particle can only be adsorbed or desorbed if at least one of its nearest neighbor sites is 
empty. This is done to mimic the short-ranged dynamical constraints in the real granular 
system. Of course, the relative magnitude e of the desorption rate and the pulse duration 
to are the parameters characterizing the process. 

Then, tapping processes have been modeled in our lattice system in the following way. 
We started from a purely random configuration, i. e. the equilibrium configuration for e = 
1/2. Then, the system was allowed to relax with e = until reaching a steady metastable 
configuration, characterized by all the holes being isolated, from which the system can not 
evolve any more. This is a convenient initial state for the compaction experiment and 
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corresponds to the loosely packed conditions used in real laboratory experiments 0. In this 
way the average initial density of particles in our tapping process has been p ~ 0.7. 

Pulses are modeled by suddenly increasing the value of e to a value greater than zero. 
This is equivalent to increase the temperature of the system. The duration of each pulse was 
to <C 1. Between pulses the system relaxes with no external excitation, i. e. with e = 0. The 
waiting time between consecutive pulses was much larger than the relaxation time needed for 
the system to become trapped in a new metastable configuration. The density was measured 
just before starting a new pulse. The whole process was designed to mimic what is done in 
real experiments. 



III. EVOLUTION WITH CONSTANT TRANSITION RATES 



In this Section we will study the evolution equations that determine the relaxation of 
the density. First, we will derive these equations for an arbitrary value of the parameter e 
characterizing the relative probability of a desorption event. Secondly, we will analyze the 
free relaxation without desorption, i. e. in the limit e = 0, and the effect of pulses separately, 
taking into account that the final state for one of the processes gives the initial condition 
for the other. 

In the following we will restrict ourselves to homogeneous and isotropic states. This 
requires to consider appropriate initial and boundary conditions, and it is consistent with 
the qualitative picture depicted in the previous Section. It will be assumed that the limit 
N — > oo has been taken. Let us define probability distributions of groups of r + 1 consecutive 
holes by 

D r (t) = (mim i+1 ■ ■■m i+r ) t = ^m t m 1+1 • ■ • m i+r p(m., t). (12) 

m 

The homogeneity of the system implies that the above expression does not depend on the 
starting site i considered. Evolution equations for the moments D r {t) are easily obtained 
from the master equation, 

^-D (t)=eD (t)-D 1 (t), (13) 



^D x {t) = e [D (t) + C ,o(*)] - D 1 (t) - D 2 (t), 



(14) 



and 



^-D r (t) = -rD r (t) - D r+1 (t) + eC , r _i(t) + e 
at 



r-l 



J'=l 



(15) 



for r > 2. Here we have introduced the probability distributions of two groups of holes 
separated by a site 



Cr,s{t) — ( m i m i+l ' ' ' m i+r rr l"i+r+2 ' ' ' m i+r+s+2)t- 



(16) 



6 



Besides, from now on we use the dimensionless time scale defined by t* = ut, although the 
asterisk is omitted for the sake of simplicity. Again as a consequence of homogeneity, the 
functions C rjS (t) do not depend on the site % taken as the origin to measure them. Moreover, 
isotropy implies the symmetry property C r)S (t) = C s>r (t). In Eq. (|T3| ) it is seen that the time 
evolution of the density of holes Do(t) involves the nearest neighbor pair distribution of 
holes Di(t). When the equation fllil) for this latter distribution is considered, the situation 
becomes more complex. Not only the three consecutive hole distribution D%(t) shows up, 
but also the second neighbor pair moment C 0)0 (t) appears. 

On the other hand, the whole hierarchy of equations gets much simpler in the limit e — > 0. 
As discussed in Sec. || this is the no desorption limit and corresponds to T — ► + (very low 
temperatures). For e = 0, Eqs. dl3|)-(|l5D reduce to 

lD?\t) = -rD?\t)-D® 1 (t), (17) 

for all r. Hereafter, the superindex indicates that a quantity is evaluated in a system 
evolving with e = 0. The hierarchy fll7|) can be easily solved by using, for instance, the 
generating function method ||10| . We introduce a generating function 

oo r 

G {0 \y,t) = j: y -Dl°\t). (18) 

r=0 r - 

From Eqs. fllTD it is obtained that G^°'(y,t) obeys the equation 

^ \y,t) + (y + l)^-G^(y,t) = 0, (19) 



whose solution is 



G<°>(^) = G [(y + l)e-*-l], (20) 
where Go(y) = G^ (y, 0) is the initial condition, that will be determined by the final situation 



after a pulse. This expression has been previously obtained in a different context fll2| . For 
large times G^ (y, t) approaches the limit 

G< >(y,oo) = G (-l), (21) 

and, consequently, 



hmDS 0) (t) = G (-l), HmDW(t) = 0, (22) 



r > 1. The last result reflects the property that for e = all the holes are isolated in the long 
time limit, i. e. they are always between two particles. Therefore, the probability of finding 
two consecutive sites with = 1 is null. This is a general property that does not depend 
on the initial conditions. Of course, the asymptotic value of = (m^}^ is determined 
by the initial state of the system, being smaller than its initial value. It must be noticed 
that the hierarchy (|T7|) admits as a stationary solution any constant value for Dq as long as 
D rst = for r > 1. 



7 



An interesting particular case is when the system is at equilibrium with a given value of 
e > before being suddenly changed to e = 0. In terms of the temperature introduced in 
Sec. [Tl]this is equivalent to a quench of the system to T = + . The initial condition for this 
process is now (see Eq. (||)), 



Then, 



and Eq. (^) yields 



D r (0) = e r+1 . (23) 
G (y) = ee^ } (24) 



(oo) = ee~ e . (25) 



For a purely random initial distribution (e = 1/2) it is Dq°^(oo) ~ 0.3033, i. e. less than 
one third of the sites are empty in the final metastable state, characterized by a "frozen" 
configuration. 

In the above discussion there was no need for considering the time evolution of the 
distributions C rjS (t) defined in Eq. ([161) . Nevertheless, it is evident that in the limit t — > oo, 
Coo approaches a constant value fixed by the initial conditions of the relaxation process, 
while CjP)(t) -> for r > or s > 0, since the last ones involve adjacent sites. 

Next we analyze the evolution of the system with e > but for a time interval to 1. 
This corresponds to the pulse preceding, and also following, each of the free relaxations 
without desorption. Therefore, the initial conditions we will be interested in correspond to 
a final state obtained after a long time relaxation with e = 0, 

£>o(0) = mo, A-(O) = for r > 1, (26) 

C ,o(0) = c , C r , s (0) = for r > 1 or s > 1. (27) 
For times t < to 1 we approximate by means of a first-order Taylor expansion using Eqs. 



D (t) ~ m + em t, (28) 

and similarly 

D x {t) ~ e (m + c )t, (29) 



D 2 (t) ~ ecot, (30) 

while distributions D r (t) with r > 3 are at least of order t 2 . Now, we can define the 
generating function corresponding to the pulse in a similar way as it was done for e = in 
Eq. (|i~8l). For this short time limit it is 
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IV. TAPPING PROCESSES 



In this Section we will use the previous results to investigate the dependence of the 
density on the number of taps. Let us consider the free relaxation with e = after the n + 1 
pulse. The initial condition for this process will be the final state reached during the pulse, 
i. e. using Eq. (|3"ID , 

+ \t^ n y 2 + O{t 2 ), (32) 

where m ,n and c 0i „ are the values of (m;) and (mjmj +2 ) at the end of the relaxation following 
the n-th tap, respectively. The time evolution of the system during the relaxation is described 
by Eq. (|20|) and in the long time limit by Eq. (pi]), that particularized for the above initial 
condition yields 

Gn+i(y, oo) - m , n - ec ,„^ (33) 

and, therefore, 

m ,„+i ~ m , n - ^et c ,n- (34) 

The above equation is expected to hold for small t but arbitrary "amplitude" of the pulses e. 
Since co, n is by definition positive it follows that the density of holes decreases and the system 
compacts monotonically as a function of the number of taps n. We stress that the density 
is measured at the end of each free relaxation as it is actually done in real experiments. 

We have checked Eq. (|34]) by comparing it with the results obtained from Monte Carlo 
simulation of the Markov process defining the dynamics of the system. An example is 
given in Fig. [1] where we have plotted both 2(p„ + i — p„)/eto and co, n as functions of the 
number of taps n. Here p n = 1 — mo in is the density of particles after the n-th tap. In 
fact, since m 0i n+i — ^o,n is a rapidly fluctuating quantity, each of the points we have plotted 
corresponds to the average of those functions over 10 consecutive taps. The data shown have 
been obtained in a system of 10 4 sites with e = 0.5 and to — 0.02, and have been averaged 
over 10 3 runs. It is seen that the prediction of the theory is verified quite accurately. For 
the sake of clarity, we have restricted ourselves to 5 x 10 3 taps, although the same behavior 
is observed until the system comes near the steady state discussed in the next section. 

Equation ([34]) indicates that the compaction process depends on the product et and, 
in that sense, eto plays in our model the same role as T in real experiments. The latter is 
defined as the ratio of the peak acceleration of the tap to the gravitational acceleration [§. 
Nevertheless, Eq. (|34|) suggests a stronger prediction, namely that the relevant time scale 
for the compaction process is r n = et^n. Of course, this will be true only if the dependence 
of co, n on n also takes place through r n , but it is easily seen that it is really so. The initial 
condition for each pulse is a trapped configuration that is metastable for e = 0. That means 
that the derivatives with respect to time of all moments are proportional to e at t — 0, and in 
the limit of short duration pulses the change in any moment in a pulse will be proportional 
to et . This proportionality is clearly kept by the free relaxation with e = that does not 



9 



introduce any new time scale in the problem. It is worth mentioning that the same results, 
i. e. Eq. (|34T), hold in the limit e — > 0, with eto <C 1. 

In Fig. |2| the relaxation of the particle density is shown as a function of the scaled time 
r n for different values of e and t . In all cases t <C 1 as required by the theory we have 
developed. The system has 10 4 sites. The initial state for all the compaction experiments 
was the same, namely the fully random distribution. As predicted, all points lie onto a single 
curve. Moreover, this scaled curve is well described by the four-parameter heuristic law 

Pn — Poo 7 \ 5 I"" J 

1 + B In f 1 + ^ 



with values of the parameters p^ = 1.10, 5poo = 0.40, B = 0.39, and r c = 3.37. We do 
not observe any dependence of these constants on the values of e or t Q , of course always in 
the limit t -C 1. The solid line in Fig. (Q) is the fit to Eq. (|35]). We have tried to derive 
analytically a logarithmic law similar to this empiric result, but we have not succeeded. It 
is not clear yet whether it is just a convenient fitting expression with four parameters or it 
has a more fundamental meaning, for instance associated to some peculiar dynamical events 
which are dominant in the relaxation of the density. In this context, it is important to 
realize that the law fails to describe the asymptotic behavior in the limit of a large number 
of taps and the steady value of the density that is eventually reached. In fact, the value 
of p^ reported above is clearly unphysical since it is larger than 1. We believe that this is 
a general limitation of the law ( |3~5"D and it is not restricted to the present model. We can 
substitute in Eq. (|35|) 

- = -, (36) 

with n c = r c /et . In this way the standard inverse logarithmic law with time measured in 
number of taps is recovered ||. But now we have an explicit dependence of n c on eto- A 
similar result was found numerically in Ref. |]] for a two-dimensional model with geometrical 
frustration. Here the dependence appears as a consequence of the relevant scale defining the 
time evolution of the system. This scale has been identified by using analytical methods. The 
value n c can be understood as the minimum number of taps needed to observe a significant 
compaction process. For n -C n c the density remains practically with its initial value. 



V. STEADY STATES 

Another point we have investigated, prompted and stimulated by the results found in 
previous works by different authors, is the possible existence of a long time steady state 
density, determined by the tapping process (i. e. the amplitude and duration of the pulses 
in the present model) but independent of the initial conditions ||. Then, we have carried 
out a series of computer experiments corresponding to the same values of to an d e but to 
different initial conditions and, in particular, to different values of the initial density. The 
states we have chosen are equilibrium states, corresponding to initial densities of particles 
0.5, 0.75, 0.9 and 0.995. The results are presented in Fig. |3|. All densities tend in the long 
time limit to the same steady value. The data in the figure have been obtained in a system 
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with eto = 0.03, but the same qualitative behavior has been found in all the studied cases. 
One important point to remark is that it is possible to start from a density higher than the 
asymptotic one and then the density decreases as the number of tapings increases. Once 
again, this behavior is analogous to what is observed in real experiments. Average densities 
above the so-called random close packing limit, which is much smaller than the crystalline 
value, are not obtained even after extensive vibratory settling. 

If we look for steady solutions of Eq. fl3H|), the only consequence we can reach is that such 
state requires Cq to be much smaller than m . A more specific statement can be obtained 
by considering the next order in the expansion of powers of to- Besides, we have simplified 
the calculations by considering a pair approximation for all the correlation functions. More 
specifically, we neglected all correlations involving more than two sites and approximated in 
Eqs. dT|) and ($L§ 



D 2 (t) = (mimi +1 m i+ 2)t 



DM 
DM 



(37) 



C ,o(t) = (mim i+ 2)t = {mim i+1 m i+2 )t + (rrii(l - m i+1 )m l+2 )t 



D\{t) [D (t) - £>!(*)]* 
D (t) 1 - D (t) 



(38) 



The above approximations can be shown to be equivalent to the dynamical mean field of 



clusters introduced by Dickman [pf] . When Eqs. (|37| ) and (^) are substituted into Eqs. 
(|13|) and (|14|), the latter become a close pair of nonlinear first order differential equations, 
namely 



dP (t) 
dt 



-D^ + eDait), 



(39) 



dt 



D (t) 



D (t) 



1 - D (t) 



(40) 



In the n+1 pulse we have to solve the above equations, for the time to 1 that the vibration 
lasts. Then, a perturbative solution in powers of to, with the initial conditions given by 



D (t = 0) = m ,„, Di(t = 0) = 0, 



(41) 



is easily obtained, 



+2 / am? 
D (t ) = m , n + et m ,n + tt ~ em o,n + e 2 «2 ,„ - 

2 V 1 - m 0,n, 



+ o(t 3 ), 



(42) 



-Di(to) = t em , n + 



em 



0,n 



m , r 



+ 



-em , n + e m ,n - 



em, 



O.ii 



em 



Q.n 



1 - m , n (1 - m , r 



+ 0(4). 
(43) 
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Afterwards, the system evolves freely with e = 0, and we measure the density of holes m ,n+i 
at the end of this relaxation. We have to solve Eqs. ( |39|) - (|40|) with e = 0. Writing them as 



(44) 



a closed second order equation for the density of holes , it is found that 



d_ 

dt 



■"-^W+inflfW 

at 



and then 



since the initial conditions for the relaxation process are the final state of the pulse, given 



by Eqs. ( fj"2]) and (|43). The value of the density at the end of the relaxation is the long time 



limit solution of the above equation, i. e. 



m , n+1 = D { °\oo) = D (t )e TO, (46) 



and using Eqs. (H2j) and 



m 0,n , 1 f . ,2 1 + m ln 



m , n +i = mo,n - tt - - h -(et ) m 0n — ^ + O(t Q m 0n ). (47) 

1 - m 0in 2 (1 - m^ n y 

As long as mo, n is much larger than eto the third term on the right hand side is negligible 
as compared with the second one, and the density of holes decreases monotonically, i. e. the 
system is compacting. To the order or approximation considered in Eq. ([|7|) a steady value 
of the density will be reached when 

m ,n+i - m 0i „ = O(^m , n )- (48) 

To get an expression for this steady density that we will denote by let assume that it 
is of the order of £q, = ct$. Therefore, 

mo,n+i - m , n = -rf^c 2 + \^l +P c + 0(t^) + O(t^) + O(t^), (49) 

and a simple dominant balance of the first and second terms on the right hand side of Eq. 
(|7|) yields (3 = 1 and c = e/2, i. e. 

4 S) = \et . (50) 

We have discarded a solution = that is always a trivial fix point for the evolution 
of the system, corresponding to all sites being occupied by particles. Let us notice that 
the above expression for is a steady solution of the evolution equations in the pair 
approximation up to and including order £g. 

In Fig. |3| we have indicated the prediction of Eq. fl5H|). An excellent agreement is found 
with the steady value reached in the simulations. A comparison for several values of eto is 
given in Fig. |j. Again the theory describes fairly well, both qualitatively and quantitatively, 
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the results of the numerical simulations. The discrepancies increases as the value of eto 
increases, as expected. This confirms that in the compaction experiment we are considering 
the steady correlations are determined mainly by the two nearest neighbor (pair) correla- 
tions. This is due to the fact that in the free relaxation processes the dynamics of the system 
is governed by the restriction that there can not be a hole next to another one, but they are 
isolated, i. e. surrounded by particles. 



VI. FINAL REMARKS 

In the framework of a simple one-dimensional lattice model with facilitated dynamics, 
we have studied the nonequilibrium evolution of a system submitted to a tapping process. 
Trying to mimic what is done in real experiments, the evolution of the system was modelled 
by a series of two alternating steps. In the first one, the system evolves perturbed by 
an external energy source, while in the second one it freely relaxes towards a metastable 
configuration. The model has been shown to share many of the characteristic features 
of granular materials under tapping. In particular, the evolution of the density can be 
accurately described by means of an inverse logarithmic law with the tapping number. 

In our model , compaction is due to the decrease of the number of holes as the system is 
being tapped. A qualitative picture of this very slow relaxation follows from the existence 
of "entropic barriers" in the system. As the number of holes lowers, the number of states 
allowing the system to relax become very small as compared with the total number of 
available states. Therefore, it is very difficult for the system to find the way to these 
bottlenecks in configuration space, being most of the time trapped exploring metastable 
configurations with almost the same density. This picture supports a strong relationship 
between structural glassy relaxation and compaction flTJJ]. In both cases there is a fast 



increase of the relaxation time of the system, becoming very large on the time scale of the 
experiment. 

Due to the simplicity of the model we have been able to obtain some detailed analytical 
results. The process is characterized by the product et , that identifies the relevant time scale 
for densification, at least in the limit of very short taps. Over this scale the time evolution of 
the system is described by a universal law, which is independent of the particular values of 
the parameters defining the system. This prediction has been compared with the numerical 
solutions obtained by Monte Carlo simulation and a very good agreement has been found. 

A main result of this paper is an analytical expression for the asymptotic density obtained 
in a tapping process. Quite surprisingly, in the limit considered there is a very simple 
proportionality relation between this density and both the duration of the taps and their 
strength, the latter being measured by the relative probability of a desorption event during 
a tap. Also this theoretical prediction has been confirmed by the numerical solution. 
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FIGURES 

FIG. 1. Plot of both co, n (solid line) and 2(p n+ i — p n )/eto (diamonds) as functions of the 
number of taps n, for e = 0.5 and to = 0.02. 

FIG. 2. Time evolution of the density of particles. Time is measured in the reduced scale 
defined in the text. In the three curves shown, the values of the parameters are e = 0.5 in all of 
them and to = 2 x 10~ 3 (diamonds), to = 0.01 (squares), and to = 0.02 (pluses). 

FIG. 3. Evolution of the density as a function of the number of taps, for four different values 
of the initial densities. The parameters characterizing the tapping process are e = 0.5 and to = 0.06 
in all cases. 



FIG. 4. Stationary value of the density of particles, p 



(») 



rn, 



as a function of the 



parameter eto characterizing the tapping process. The diamonds are numerical results, while the 
solid line is the analytical expression given by Eq. (50). 
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